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Abstract 

We study the classical dimer model on a square lattice with a single vacancy 
by developing a graph-theoretic classification of the set of all configurations 
which extends the spanning tree formulation of close-packed dimers. With 
this formalism, we can address the question of the possible motion of the 
vacancy induced by dimer slidings. We find a probability 57/4 — lOv^ for 
the vacancy to be strictly jammed in an infinite system. More generally, the 
size distribution of the domain accessible to the vacancy is characterized by a 
power law decay with exponent 9/8. On a finite system, the probability that 
a vacancy in the bulk can reach the boundary falls off as a power law of the 
system size with exponent 1/4. The resultant weak localization of vacancies 
still allows for unbounded diffusion, characterized by a diffusion exponent that 
we relate to that of diffusion on spanning trees. We also implement numerical 
simulations of the model with both free and periodic boundary conditions. 



1. Introduction 

1.1. The problem 

The statistical mechanics of rigid dimers arranged on a lattice is relevant to many 
physical systems and has appeared repeatedly in the literature over the years. The study 
of dimer models can shed light on diatomic gases, the thermodynamics of adsorbed films 
1^ and the classical limit of resonating valence bond (RVB) models of high-temperature 
superconductivity 0. 

The case of close-packed dimers on planar lattices is exactly solvable and the associ- 
ated mathematical techniques have proven very powerful [^If. It provides a paradigm 
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Fig. 1: An example of a dimer move on the square lattice. The dimer 
slides into the empty site (vacancy) of configuration (a), resulting in a new 
configuration (b) where the vacancy has jumped by two lattice spacings. 

of a geometrically constrained statistical system with a strong interplay between the 
physical degrees of freedom and the symmetry of the underlying lattice. 

New challenges appear if we extend the close-packed dimer models by allowing 
for defects in the form of vacancies, namely sites not covered by dimers ||7|,p|-p!T| . The 
presence of defects allows dynamical moves for dimers by sliding them into empty sites, 
as illustrated in Fig. 1. This in turn induces motion of the vacancies in the form of 
discrete jumps, each by two lattice spacings. It is of great interest to characterize as 
explicitly as possible the nature of the constrained dynamics of vacancies as a model 
glassy system. In particular, we would like to determine the extent to which a vacancy 
can diffuse, both spatially and temporally. 

In this paper, we consider the simplest case of an isolated vacancy in a sea of dimers 
on the square lattice. Our primary interest is the structure of the space accessible to the 
vacancy. In particular, we address the question of whether or not the vacancy is local- 
ized, namely confined to a finite region. We find that our model exhibits localization but 
in a very weak form. Although the motion of the vacancy in a fixed dimer background 
is localized to a finite domain, the mean size of this domain nevertheless diverges upon 
averaging over all possible dimer backgrounds. We call this property weak localization, 
as it allows vacancies to diffuse arbitrarily far on average. It is also important to in- 
vestigate the kinetics of vacancy diffusion. We find an anomalous diffusion exponent 
determined both by the internal structure of the space available to the vacancy and the 
effects of weak localization which control the size of this space. 



1.2. Outline and summary of results 

Our paper is structured as follows. In Sect. 2, we recast our dimer problem with 
a single vacancy as a model of spanning graphs on a rectangular grid by extending 
a famous construction by Temper ley. This construction is recalled in Sect. 2.1 and 
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provides a bijection between dimer configurations with the vacancy on the boundary 
of the grid and spanning trees on a subgrid. The generahzation of this construction, 
in Sect. 2.2, to a vacancy in the bulk leads to more general spanning webs, consisting 
of a central tree component on which the vacancy can freely diffuse, surrounded by a 
number of nested loops with branches, which act as cages for the vacancy. We then 
derive in Sect. 2.3 a determinant formula for the partition function of these spanning 
webs. 

Sect. 3 is devoted to a finite size study of dimers on a square grid of size L, with a 
vacancy at the center, and the asymptotic limit L —> oo. In Sect. 3.1, we use our deter- 
minant formula to compute the probability that the vacancy can reach the boundary of 
the grid, in which case it can reach all sites of the grid. This "delocalization probability" 
is found to decay with the system size as L~^/'^. In Sect. 3.2, we analyze the distribution 
p{s) for the size s of the domain accessible to the vacancy in an infinite system. In the 
spanning web formulation, this domain is nothing but the central tree component that 
contains the vacancy. We characterize p{s) by first giving the exact value of which 
measures the probability that the vacancy is strictly jammed. We then derive its large s 
behavior p{s) ~ from the associated scaling of the delocalization probability. This 

power law behavior with an exponent larger than —2 is responsible for the announced 
weak localization property. We finally discuss in Sect. 3.3 the diffusion exponent on 
spanning webs (ry) and on spanning trees (ryo) and find the relation ry/ryo = 7/8, which 
measures the slowing of diffusion by weak localization. 

We complete our analysis in Sect. 4 by a number of numerical simulations both 
on the statistics of spanning webs and the dynamics of diffusion. In Sect. 4.1, we show 
how to modify the well known Propp- Wilson algorithm for the generation of random 
spanning trees so as to obtain spanning webs with a uniform measure. We present 
results for the delocalization probability and for the size distribution |>(s), finding good 
agreement with the analytic predictions. We check the universality of p{s) by also 
carrying out simulations with periodic boundary conditions. There we used the efficient 
"pivot algorithm" to generate configurations directly in the dimer formulation. Sect. 4.2 
is devoted to the numerical computation of the diffusion exponents ryo on spanning trees 
and ?y on spanning webs, whose ratio agrees with the analytic prediction. 

Appendix A gives a heuristic derivation of the decay exponent for the delocalization 
probability via Coulomb gas arguments. Appendix B discusses the technical details of 
the extension of the Propp- Wilson algorithm to the generation of spanning webs. 



2. Vacancy in a dense sea of dimers: graph-theoretic treatment 
2.1. Tree formulation 

Let us first recall Temperley's bijection between dimer configurations on a rectan- 



gular grid and spanning trees on a grid with double mesh size [|T2| , p!3[] . More precisely, 
let us consider a rectangular grid with 2L + 1 columns and 2M + 1 rows, so that the 
total number of vertices is odd. The vertices of the grid can be colored, say in black 
and white, so that neighboring vertices have different colors, with the vertex in the 
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Fig. 2: A 9 X 7 grid (a) is bicolored with one more white than black vertex. 
The set of white vertices is made of two intercalated grids of double mesh 
size (b) with respective sizes 5x4 and 4x3. The first one (c) will be 
referred to as the odd white grid and the second one (d) as the even white 
grid. The latter is extended (dashed lines) so as to include an extra vertex 
(here represented by a double boundary) dual to the exterior face. 



upper right corner being white. With this coloring, there is one more white vertex than 
black vertex and the grid can be fully covered by a set of dimers with a single vacancy 
on a white vertex. Note that the set of white vertices can be viewed as made of two 
intercalated grids, of size (L + 1) x (M + 1) and L x M respectively, both with double 
mesh size (see Fig. 2). We shall refer to these grids as the odd white grid and the even 
white grid respectively. It is convenient to extend the even white grid by an additional 
white vertex dual to the exterior face of the original grid, together with edges from that 
new vertex to all vertices on the boundary of the even white grid (see Fig. 2). With this 
convention, the "extended" even white grid is simply the dual graph of the odd white 
grid. Note finally that the black vertices sit precisely on the edges of either of these two 
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Fig. 3: A sample dimer configuration on a 17 x 13 lattice (a) with a single 
vacancy ((^) in the upper right corner. In (b) we keep only those dimers 
that cover a vertex of the odd white grid and translate them in (c) into 
oriented edges of that grid, resulting in a spanning tree configuration whose 
edges are oriented toward the upper right corner. In (d), we also show the 
dual tree (red) which spans the extended even white grid and is oriented 
to the exterior face. The configuration (a) is recovered by replacing every 
oriented edge of (d) by a dimer on its first half. 

dual white grids. 

Let us assume that we place the vacancy on a white vertex of the boundary of the 
grid, necessarily part of the odd white grid. Now any other vertex of the odd white grid 
carries a dimer (see Fig. 3-(a)). This dimer selects an edge of the odd white grid, which 
we decide to mark and orient from the white vertex carrying that dimer to its neighbor 
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on the odd white grid (see Fig. 3-(b,c)). The graph made of aU these marked oriented 
edges together with all the original vertices of the odd white grid forms a spanning tree 
of the odd white grid whose edges are moreover oriented toward the vacancy (see Fig. 3- 
(d)). To understand the absence of loops in the graph, note that a loop would enclose 
an interior region of the original lattice with an odd number of vertices that therefore 
could not be fully covered by dimers (here it is crucial that we have taken the vacancy 
to lie strictly on the boundary). As the graph has exactly one more vertex than edges 
(as the vertex carrying the vacancy does not give rise to a marked edge) , it is necessarily 
made of a single tree spanning the whole odd white grid. The edge orientations produce 
the unique flow on the tree toward the vacancy (with exactly one edge exiting from all 
vertices but the vacancy). In particular, if we fix the position of the vacancy, say in 
the upper right corner, the above construction provides a bijection between fully-packed 
dimer configurations of a (2L + 1) x (2M + 1) grid with one vacancy in the corner 
and spanning trees of the associated (L + 1) x (M + 1) odd white grid. To go from the 
spanning tree configuration back to the original dimer configuration, we first orient each 
edge of the tree so as to reproduce the unique flow toward the upper right corner and 
replace each edge of the tree by a dimer on the first half of the edge. We then consider 
the "dual tree" of the spanning tree, made of those edges of the extended even white 
grid that do not cross the edges of the spanning tree. This dual graph is itself a tree 
spanning the extended even white grid which we orient toward the exterior vertex (see 
Fig. 3-(d)). We finally repeat the above construction and replace each oriented edge of 
the dual tree by a dimer on the first half of the edge. 

Clearly, we have a similar bijection if we place the vacancy at some arbitrary 
but fixed white vertex on the boundary of the grid. This simply amounts to consider 
spanning trees with another orientation of the edges, now pointing toward the new 
position of the vacancy on the boundary. If we leave the position of the vacancy free 
along the boundary, we clearly get a 2(L + M) to 1 mapping instead. 



Fig. 4: The elementary motion of the vacancy resulting from the sliding 
of a dimer (top). In the tree formulation, this amounts to reversing the 
orientation of the corresponding edge on the (oriented) spanning tree. 

Let us now consider the motion of the vacancy generated by sliding of dimers. In the 
tree formulation, as illustrated in Fig. 4, performing an elementary slide corresponds to 
picking an edge adjacent to the vacancy (and necessarily pointing to it) and reversing its 
orientation so that the fiow now points to a new vertex. We can then repeat the process, 
which allows the vacancy to reach any vertex on the spanning tree, hence any vertex 
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of the odd white grid. In the tree formulation, the dimer configurations accessible by 
slidings arc therefore described by the same spanning tree and differ only by the choice 
of the (arbitrary) vacancy site on the odd white grid toward which we orient all the 
edges of this tree. Each orientation may then be mapped into dimers as above. Note 
that the dimers corresponding to edges of the dual tree cannot be affected by slidings 
and are in practice frozen. Note finally that the property that the vacancy can reach 
any site of the odd white grid holds only because the vacancy was on the boundary in 
the first place. As shall see just below, this property is not true in general for vacancies 
that lie inside the grid. To conclude, the motion of the vacancy can be analyzed as 
simple diflFusion on a spanning tree picked uniformly at random. 

2.2. Extension to webs 

Let us now consider the more general case where the vacancy originally lies on an 
arbitrary vertex of the odd white grid, not necessarily on the boinidary (see Fig. 5- (a)). 
We can repeat the above construction by marking and orienting, for each vertex of 
the odd white grid except the vacancy vertex, the edge that is selected by the dimer 
covering that vertex (see Fig. 5-(b)). As before, we consider the graph made of all these 
marked oriented edges together with all the original vertices of the odd white grid. As 
there is exactly one oriented edge exiting from each vertex but the vacancy vertex, 
any connected component of this graph is either a tree containing the vacancy vertex 
and with edges pointing to that vertex or it is made of an oriented loop with attached 
branches oriented toward it (see Fig. 5-(c)). In this latter case, the connected component 
cannot contain the vacancy vertex. Moreover, as before, any such loop encloses an odd 
number of vertices, so it must encircle the vacancy., which is possible only if the vacancy 
does not lie on the boundary. To summarize, the general structure of the graph is a tree 
containing the vacancy vertex (the tree may be as small as a single vertex) surrounded 
by a number of nested loop components (made of a loop and attached branches) so that 
the graph spans the entire odd white grid. We shall call such a structure a "spanning 
web" of the odd white grid rooted at the vacancy vertex. 

Conversely, given such a graph, we note that the "dual web" made of those edges 
of the extended even white grid that do not cross the edges of the spanning web is itself 
a spanning web of the extended even white grid (see Fig. 5-(d)). It is made of a set of 
nested loop components, with branches, complemented by a tree rooted at the exterior 
vertex. The loops of the dual web and that of the original spanning web alternate 
and therefore match in number. We can translate the spanning web configuration into 
a dimer configuration by successively orienting the edges of its tree component toward 
one of its vertices, where we put the vacancy, the edges of each loop in the same direction, 
chosen arbitrarily for each loop, and finally the edges of the attached branches toward 
the loops. Similarly, we orient the edges of the tree component of the dual web toward 
the exterior vertex, those of the loops in the same direction, arbitrarily for each loop and 
those of the attached branches toward the loops. The dimer configuration is obtained 
by replacing each of these oriented edges by a dimer on the first half of the edge. 

The above construction provides a bijection between dimer configurations with 
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Fig. 5: A sample dimer configuration (a) with a single vacancy ((^) in the 
bulk of the odd white grid. In (b) we show only those dimers that cover 
a vertex of the odd white grid. When mapping them into oriented edges 
of that grid, the resulting spanning web (c) consists of a tree component 
containing the vacancy, surrounded by a number of loop components (here 
2). The edges belonging to branches of the tree component are oriented 
toward the vacancy while those belonging to branches of the loop compo- 
nents are oriented toward the loop. On each loop, all edges have the same 
orientation. In (d), we also show the "dual web" (red) which spans the 
extended even white grid and is made of an equal number of loop com- 
ponents surrounding the vacancy together with an extra tree component 
attached to the exterior face. The edges belonging to branches of this tree 
component are oriented toward the exterior face while those belonging to 
branches of the loop components are oriented toward the loop. Again, all 
edges have the same orientation on each loop. The configuration (a) is 
recovered by replacing every oriented edge of (d) by a dimer on its first 
half. ^ 



a vacancy at a fixed position on the odd white grid and spanning webs whose tree 
component contains the vacancy, together with a choice of orientation for each loop of 
the spanning web and each loop of the dual web. As there are 2 possible orientations 
per loop and an equal number of loops on both web configurations, we can get rid of 
the orientations by counting each loop of the original spanning web with a degeneracy 
factor 4. 

Again, let us examine how the vacancy can move under the sliding of dimers. As 
before, an elementary slide amounts in the spanning web language to reversing the 
orientation of one edge pointing to the vacancy so that it points to a new vertex. Under 
repeated moves, the vacancy can reach every white site of the tree component of the 
spanning web, which therefore constitutes a complete specification of the set of sites 
accessible to the vacancy. Note that all the dimers corresponding to oriented edges of 
either the loop components of the spanning web or any component of its dual web are 
frozen. In particular, the vacancy cannot cross a loop of either web. 

On a finite grid, we may therefore divide the configurations into two classes, de- 
pending on whether or not the vacancy can reach the boundary. If the vacancy can reach 
the boundary, there cannot be any loop as the loops are required to encircle the tree 
component of the spanning web. In this case, the spanning web reduces to a spanning 
tree and the vacancy can reach any vertex of the odd white grid. If the vacancy cannot 
reach the boundary, there must be a loop component containing all the white boundary 
vertices. There are in general several nested loops, the interior-most one acting as a cage 
for the vacancy. The motion of the vacancy may then be analyzed as simple diffusion 
on the tree component of the web. For a fixed initial position of the vacancy on the 
grid, the ratio between the number of dimer configurations for which the vacancy can 
reach the boundary and the total number of dimer configurations can be interpreted as 
a delocalization probability. This probability will be studied in detail in Sect. 3.1 below 
for the case of a square grid with vacancy in the center. 

On an infinite grid, the vacancy is always localized to a finite tree component but 
we shall see that the mean size of this tree in fact diverges. 



2. 3. Determinant formulae 

It is well known that the number of spanning trees on a graph with n vertices is 
given by any principal (n — 1) x (n — 1) minor of the Laplacian matrix of this graph. 



This is the celebrated matrix-tree theorem attributed to Kirchhoff [0. Recall that the 
coefficients of the Laplacian matrix simply read 

{di for i = j 
— 1 for i and j neighbors (2.1) 
otherwise , 

where i and j are vertices of the graph and di is the degree of i (number of incident 
edges) on the graph. 

Here we are interested in the case of a rectangular grid where di can be 2, 3 or 4 
according to whether i lies at a corner, on the side or in the bulk of the grid. Computing 
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the minor amounts to removing a given vertex io, i.e. restricting the indices of the matrix 
to aU i ^ iQ. Note that although zq is removed, the degrees of its neighbors are kept 
unchanged. The number Ztree of spanning trees therefore reads 

Ztree = det(A,,)^^^.^^^ . (2.2) 
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Fig. 6: A graphical representation of a term in the expansion of the de- 
terminant in Eq. ( |2.2|) by use of oriented blue and red edges. There is 
exactly one outgoing arrow from every vertex distinct from io. Blue edges 
must form oriented cycles, each cycles receiving a weight —1. Red edges 
may either form oriented cycles or branches oriented toward a (blue or red) 
cycle or toward ig. 



Let us sketch a proof of this result, along the lines of Ref. [T3|, that can be easily 
extended to spanning webs. Writing the determinant as a sum over permutations, we see 
that the only permutations having a non-zero contribution are those with the following 
two types of cycles: (i) trivial cycles of length one (fixed points), each contributing a 
weight equal to the degree the associated vertex, (ii) cycles of length larger or equal 
to 2 for which any two successive elements of the cycle are neighbors on the grid. 
The net contribution of any such non-trivial cycle to the determinant (including the 
signature of the permutation) is easily seen to be —1. It is convenient to have a pictorial 
representation of each such permutation as follows. Any cycle of type (ii) above can 
be represented as a closed oriented loop on the grid by joining each vertex to its image 
under the permutation. The loops moreover avoid the removed vertex zq. Such loops 
are represented in blue in Fig. 6 and will be henceforth referred to as blue loops. Each 
blue loop contributes a factor —1. For the trivial cycles (i), the contribution di of the 
associated vertex i is properly accounted for by considering the di possible choices of an 
edge incident to that vertex. For each choice, we mark and orient the edge away from 
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the vertex. Such oriented edges are represented in red in Fig. 6 and wiU be referred 
to as red edges. Any permutation is thus represented by a number of configurations 
(corresponding to the different choices of red edges for trivial cycles) made of blue loops 
and red edges with the only constraint that there is exactly one outgoing (blue or red) 
edge from each vertex of the grid but the removed vertex zq- Each configuration is now 
weighted by a factor 1 per red edge and —1 per blue loop. Now clearly, any configuration 
in which the red edges form a loop is canceled exactly by a similar configuration where 
all the edges of this loop are blue. The only configurations that remain therefore consist 
of red edges forming a spanning tree, necessarily oriented toward the removed vertex 
zq. This completes the proof. 

From the determinant formula in Eq. (2^), one can easily derive a closed formula 
for Zti-ee as a product over eigenvalues of A [|TB|. 
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Fig. 7: The odd white grid (a) completed with a seam from io to the 
boundary (thick dashed line). This seam crosses a number of edges that 
we call seam edges and that we orient so as to point to the right when 
fiowing on the seam from io to the boundary. In the graphical expansion 
of the determinant (b), the seam modifies the weight of the blue edges 
crossing it, resulting in a new weight for blue loops. Here the blue loop 
that winds around io passes through exactly one seam edge (with the 
wrong orientation), resulting in a weight factor — 1/a instead of —1 in the 
absence of a seam. On the contrary, the lower-left blue loop that does not 
wind around io passes exactly twice through seam edges, with canceling 
orientations so that its weight is unchanged. 



The expression ( |2.2| ) is easily modified to obtain the number Z„eb(^o) of spanning 
webs on a rectangular grid rooted at some fixed vertex io. We simply construct an 
oriented seam from io to the exterior face that crosses edges only (see Fig. 7). Each 
crossed edge, which we denote a seam edge, can be oriented so as to point right when 
fiowing along the seam. We may then modify the Laplacian matrix into a matrix A (a) 
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defined as: 



A{a)ij = < 



di for i = j 

for i, j neighbors not separated by a seam edge 
for i, j neighbors separated by a seam edge oriented from i to j 
for i, j neighbors separated by a seam edge oriented from j to i 
otherwise . 

(2.3) 

For any choice of the seam, the quantity 



-1 

-a 
1 

a 





Z{y;io) = det{A{a)ij),^.^,^ (2.4) 
then counts the number of spanning webs with root zq and with a weight 

y = 2 - a - - (2.5) 

o 

per loop. In particular, we recover Ztree = Z{0; ig) for a = 1 while for spanning webs, 
the desired weight 4 per loop is obtained by choosing a = — 1, i.e. 



Z^eb(^o) = Z{A;io) = det (A(-l),,)^^^._,^^ . (2.6) 

To obtain this result, we simply use for the determinant the same representation as above 
in terms of blue loops and red edges. Defining the algebraic number of seam crossings 
as the number of passages across seam edges in the correct orientation minus that in 
the wrong orientation, this number is zero for a loop that does not encircle io, +1 for 
a loop that encircles zq clockwise, and —1 for a loop that encircles counterclockwise. 
The blue loops that wind around the vertex io now get a modified weight —a or — 1/a 
according to whether they are oriented clockwise or counterclockwise. As the weight of 
red edges is unaffected, any clockwise (resp. counterclockwise) red loop encircling the 
vertex zq combined with the same configuration with a blue loop instead gives rise to 
a weight 1 — a (resp. 1 — 1/a). Summing over both orientation results in the weight y 
above. 



3. Finite size analysis 

3.1. Determinant calculations and asymptotic estimates 

In this section, we shall restrict our analysis to the case of dimers on a square grid 
with the vacancy in the center. We shall consider only squares of size (4L + 1) x (4L + 1) 
so that the center vertex zq lies on the odd white grid, of size (2L+1) x (2L+1). With this 
geometry, we can rely on the fourfold rotational symmetry of the problem to perform a 
block-diagonalization of A (a) which translates into the factorization: 

Z{y; io) = Pl {a)PL {-a)PL {ia)PL (-i a) , (3.1) 
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where = a with a as in Eq. (|2.5| ), namely 

y = 2-a^-^. (3.2) 

Here PL(a) is the determinant of a square matrix A(a) of size L(L + 1) (compared to 
4L{L + 1) for A(a)) corresponding to a modified Laplacian on a quarter grid. More 
precisely, we consider the graph of Fig. 8 obtained by keeping only one quadrant of the 
original odd white grid and completing it by adding L oriented "winding edges" joining 
the vertices of two consecutive sides of the quadrant as shown. The matrix A(a) has 
diagonal elements A{a)ii = di where di is the degree of site i on the graph of Fig. 8 
completed by the vertex iq as shown, which is also the degree of the corresponding 
vertex on the original odd white grid. The off-diagonal elements K{a)ij are given by 
the sum over all edges connecting i and j of a contribution equal to —1 for regular 
(non- winding) edges, —a (resp. — 1/a) for winding edges oriented from i to j (resp. 
from i to i). 

It is interesting to note that Phiot) can be interpreted as the number of spanning 
web configurations rooted at the center of the square grid of size (2L + 1) x (2L + 1) 
and which are symmetric under 7r/2 rotations around the vacancy vertex. For any 
such symmetric spanning web, the loops now come with a weight y = 2 — a — 
In particular, for a = —1, Pl(— 1) counts the number of fourfold symmetric dimer 
configurations on the square grid of size (4L + 1) x (4L + 1). 

Using Mathematica, we have computed Pl(ck) exactly up to L = 60, corresponding 
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Fig. 9: The numerical value (a) of log]^o('^i-,i) versus z for L = 10, 20, 30, 
40, 50, and 60, from bottom to top. All the data fall on the same scaling 
curve (b) upon using reduced variables i/L and log]^Q(7ri^^i)/ log]^Q(7ri;,^o)- 



to dimers on a 241 x 241 grid. For L = 1, 2, 3, we have 



Pi (a) = 4 - a + - 
a 



n / 2 1 



P2{a) = 178 - 60 a + - + + ^ 



a 



P3(a) = 82128 - 31667 a + - + 1160 + 



1 



a 



a" + 



(3.3) 



More generally, PL(ct) is a Laurent polynomial of a. of the form 



-a] 



i=-L 



(3.4) 



where the til,! are positive integers and satisfy TiL,-i = '^L,i- Their numerical values for 
L = 10, 20, 30, 40, 50 and 60 are displayed in Fig'. 9. 

From our data, we first conjecture the amusing fact that -PL(i) counts the number 
of fully packed dimer configurations on a cylinder of height 2L and circumference 2L + 1 . 

More relevant to our study, we expect a large L behavior 



PL(a) ~ c(a) 



(3.5) 



where the area entropy factor /i and the boundary entropy factor A are independent of 
a while the prefactor c and the exponent 7 depend on a. The values of /x and A are 
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known to be 
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H = exp (^-^ j dk j (i£ In (4 - 2 cos /c - 2 cos t)^ 

= exp = 3.209912300728158 . . . '^^'^^ 

A = 72-1 , 



where G = X]i^o(~^)V(2^ + l)^ is Catalan's constant. We have checked the consistency 
of our exact finite size data with these values. For instance, the value of /i can be es- 
timated by considering the ratio ^^(a) = Pi_|_2(a)PL(tt)/(-PL+i(«))^5 which converges 
to at large L, and applying to Al standard convergence acceleration techniques. In 
this paper, we found it convenient to use a simple linear convergence algorithm which 
consists, for any sequence Ul tending at large L to C/qo, in building new sequences 
= A^^^ ^L'^UL/kl'j for increasing integers k, where A'^'^^ denotes the k-th iteration 

of the finite difference operator A/(L) = /(L + l) — /(L). The sequence Uj^^ is expected 
to tend faster to Uoo for larger k as long as k is kept reasonably small (in practice we 
used mostly k = 4 and went as far as /c = 11 for the jamming probability of Sect. 3. 2 
below). For instance, if Ul has a Taylor expansion in 1/L, then the Taylor expansion 

of Uj^^ has its first k correction terms vanishing. 

3.211 
3.2105 
^ 3.21 

S-l 

^ 3.2095 
I 3.209 
W 3.2085 
3.208 
3.2075 

10 20 30 40 50 60 

L 




Fig. 10: Estimate for the area entropy // at a = ±1 from the data of 
for L up to 60. 

The estimates for /i when a = +1 and a = —1 are plotted in Fig. 10 as a function 
of L. We obtain the value /x = 3.2099123(1), fully consistent with the exact analytic 
expression above. 
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More interesting is the value of the exponent 7. Here we shaU be interested 
only in real values of y in the range [0, 4] . This corresponds to taking complex val- 
ues of a on the unit circle. The value of 7 can be estimated from the quantity 
Bl = -iLy2)\n {PL+3{a)iPL+i{a)f/{PL{a){PL+2{a)f)). which converges to 7(a) 
at large L. Again the results are improved by use of our convergence acceleration 
algorithm. 
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Fig. 11: Estimates of the exponent 7(0;) for, from top to bottom, a = 1, 
a = exp(±i7r/4), a = ±i, a = exp(±3i7r/4) and a = —1, as obtained from 
the data of PL(ct) for L up to 60. 



Fig. 11 displays the estimated values of 7 as a function of L for the various values 
of a corresponding to y = (ct'* = 1) and y = 4 (a^ = —1). We obtain the asymptotic 
estimates 



7(1) = 0.749999(1) 

= -0.249999(1) 
7(i) = 7(-i) = 0.000000(1) 



(3.7) 



for the four roots of = 1 and 



7(exp(i7r/4)) = 7(exp(-i7r/4)) = 0.3125002(2) 
7(exp(3i7r/4)) = 7(exp(-3i7r/4)) = -0.187499(1) 



(3.8) 
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for the four roots of = 



— 1. From these estimates, we conjecture the exact values 



7(1) =3/4 
7(-l) = -1/4 

7(i)=7(-i) = (3.9) 
7(exp(i7r/4)) = 7(exp(— i7r/4)) = 5/16 
7(exp(3i7r/4)) = 7(exp(-3i7r/4)) = -3/16 . 




0.5 1 1.5 2 



e 

Fig. 12: Estimated value of 7(a) for a = exp(i7re) with < e < 2, 
as obtained from the data of consecutive values of Pl(«) for L near 60 
(red line). The dashed green line indicates the conjectured exact result 
7(a) = 3/4-e(2-e). 

More generally, Fig. 12 displays the estimate of 7(a) for all values of a on the unit 
circle. These data are fully consistent with the following analytic expression for 7: 

3 

7(a) = - - e(2 - e) , a = exp(i7re) , 0<e<2. (3.10) 

This form for the exponent, quadratic in the phase e, is frequently encountered in models 
with statistical weights associated to loops and will be assumed exact in all subsequent 
analysis. 

From Eq. ( |3.1| ), we immediately deduce from our data the exact generating function 
Z^iy) = Z{y; io) for a (2L + 1) x (2L + 1) odd white grid with io in the center for L up 
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to 60. For instance, we get 
Zi(y) = 192 + y 

Z2{y) = 557568000 + 10474560 y + 

Zsiy) = 19872369301840986112 + 647704492383277056 y + 1642581444224^2 + . 

(3.11) 



Note that Zqo{0) and .^60(4) are 7323 digit numbers! From Eqs. ( p?T|) and (|3.5| ), we 
have the foUowing large L behavior for Z{y; io) 



Z{y;io) ~ c{y) 



l^i2L+lf^4L 



(3.12) 



where 7(1/) = 7(a) + 7(— a) + 7(i«) + 7(— ia) with a as in Eq. Q and c{y) 
c{a)c.{—a)c{ia)c.{—ia)/ii. 
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Fig. 13: A direct estimate of the exponents 7(0) and 7(4) obtained from 
the data of Zl{0) and Zl{4) for L up to 60. 

In particular, for y = {a^ = 1) and y = 4 (a^ = —1), we obtain 7(0) = 3/4 
1/4 + + = 1/2 and 7(4) = 5/16 + 5/16 - 3/16 - 3/16 = 1/4, namely: 



Ztree ~ c(0) 

Z^ehiio) ~ c(4) 



^(2L+1)^^4L 

Li/4 



(3.13) 



Note that the first exponent (1/2) can be obtained directly from the asymptotics of the 
closed product formula for Ztree [H^- As shown in Fig. 13, these values are corroborated 
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by a direct estimate of 7 from the exact values of ^l(O) and Zl{4) for L up to 60. We 
indeed estimate 7(0) = 0.499999(1) and 7(4) = 0.250001(1). 

Taking the ratio of Z^ree and Zweb(^o)5 we deduce that the delocaUzation probabihty 
Vl of a vacancy at the center of a (4L + 1) x (4L + 1) square decays at large L as 

- ' (3-14) 

with some prefactor c(4)/c(0). 

More generally, we immediately get from the expression ( |3.10| ) (with a restricted, 
without loss of generality, to0<e<l/2): 

7(2/) = ^^ + l : 2/ = (2 cos{nu)f , -^<u<^, (3.15) 

where we have set u = 1/2 — 2e. This expression can be alternatively obtained from a 
heuristic argument based on an effective Coulomb gas description of the model. This 
argument is presented in Appendix A. 

The degree of localization of the vacancy is also measured by the number C of loops 
of the spanning web that surround it. The average of C over all spanning webs reads 

{C)=y-^\nZ{y;io) . (3.16) 
ay 

At large L, we find from the asymptotic behavior (|3.12| ) that 



n 11 
(£) ~ -y— 7(y) InL = InL , (3.17) 

ay 7rtan(7rw) 

where we have used the explicit form ( |3.15| ) of ^{y) with y — {2 cos(7rw,))^. For instance, 
when y = 4, (i.e. u — > 0), we find that the average over all possible dimer configurations 
of the number of loops around the vacancy scales as 

(£)~^lnL (3.18) 

at large L. 

From the exact values of Zi(y), we also have access to the exact probability for 
spanning webs to have C loops on a finite square grid with the vacancy in the center. 
Fig. 14 displays these probabilities in the case y — 'i for £ = (which is nothing but Vl) 
and C = 1,2, with L up to 60. Note that for this range of grid sizes, the probability that 
£ > 3 is negligible (less than 10~^). All these probabilities are expected to eventually 
decay as L~^/^ at large L. 
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Fig. 14: The exact probabilities for having C loops in a spanning web 
with y = 4 for an odd white grid of size (2L + 1) x (2L + 1) with the 
vacancy in the center. We represent only the cases of £ = 0, 1, 2 as the 
other probabilities are negligible. 



3.2. Asymptotic size distribution 

We have established so far that, for a finite grid, the delocalization probability Vl 
decreases with the system size L as a power law with exponent 1/4. Let us now see 
how to extract from this result the value of two other exponents which characterize 
the possible motion of a single vacancy in a sea of dimers on an infinite grid. A first 
exponent, S, characterizes the size distribution of sites accessible to a vacancy. More 
precisely, let us consider again a grid of size (4L + 1) x (4L + 1) covered by dimers with a 
vacancy at the center and consider the probability that the tree component of the 

associated spanning web has s vertices. We expect this probability to tend at large L 
to some limiting law p{s) = limx^^oo Pl('S), with a finite value for all positive integers s. 
This asymptotic distribution should be universal in the sense that it should not depend 
on the precise initial position of the vacancy in the bulk or on the imposed boundary 
conditions. 

Of particular interest is the value of p{l) which measures the probability that the 
vacancy is fully jammed. This value is easily estimated from the exact values of pl(1) 
at finite L. For s = 1, the tree component of the spanning web is reduced to the single 
vertex io. In the determinant formulation of Sect. 2.3, this simply means that no red 
edge points toward io. Imposing s = 1 therefore amounts to changing di di — 1 in 
Eq. (|2.3|) for the diagonal terms associated with the four neighbors of io on the grid. 
Again we can rely on the fourfold symmetry to reduce by four the size of the matrix 
involved. We have computed Pl(1) for L up to 50, as plotted in Fig. 15. By applying 
our convergence algorithm, we estimate p(l) = 0.10786437626904951198(1) from which 
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Fig. 15: The probability that the vacancy at the center of a (4L + 

1) X (4L + 1) grid is strictly jammed. 



we conjecture the exact value 



p{l] 



57 

T 



10 



(3.19) 



This value was identified thanks to Plouffe's Inverter ||T^ applied on the first 10 digits 
of 1/ \/p{l). We then verified that it indeed reproduces all 20 digits of p{l) above ! It 
would be nice to have an analytic proof of this result. 

We have also computed Pl(2) exactly for L up to 34, from which we estimate 
p(2) = 0.055905353801942(1). We have not been able to guess an exact expression for 
p{2). 

At large s, the distribution p{s) should behave as 



p{s) 



(3.20) 



with an exponent 5 that we shall now compute. First, since the delocalization probability 
Vl tends to at large L, it follows that the vacancy in an infinite system is always 
localized. This implies that the distribution p(s) is normalized to 1, namely P{^) — 

1. 

To evaluate 5, we consider the probability 1 — Vl that the vacancy in a finite grid is 
localized, i.e. cannot reach the boundary of the grid. At large L, we can indeed estimate 
this probability as the probability that s remains less than a maximal size of order 4L^, 
namely 

4L^ oo „oo 

1 ^ ^p(s) = 1 _ p(s)ds ~ 1 -const. L^^i--^) (3.21) 
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Comparing with ( p.l4|) , one finds 



2(1 - 5) = -1, (3.22) 



from wliidi we deduce 

S=l. (3.23) 

Despite tfie fact that the vacancy is locahzed on a finite tree, the result that 6 < 2 
implies that the mean size (s) of this tree nevertheless diverges. For some observables, 
the vacancy behaves on average as if it were delocalized. In this sense, our model 
possesses a rather unusual property which we may call weak localization. 

3.3. Diffusion exponent 

As an illustration of this unusual property, let us now discuss the dynamics of 
diffusion of a single vacancy in a dense sea of dimers. As we discussed, the motion of 
the vacancy is induced by a dimer sliding axially into it, resulting in a jump of the 
vacancy by two lattice spacings. In the spanning web language, this simply amounts to 
moving the vacancy to one of its neighbors on the tree component of the spanning web. 
A natural choice of dynamics consists in choosing, at each time step, one of the four 
neighbors of the vacancy on the odd white grid at random, and moving the vacancy 
to that neighbor if possible, i.e. if this neighbor belongs to the tree component of the 
spanning web. We then define S{t) as the total number of sites visited by the vacancy 
after t steps. Although S{t) is clearly bounded by the size s of the tree component of 
the spanning web, we expect from our weak localization property that it can become 
on average arbitrarily large, with a large time behavior of the form 

{S{t)) ~ W , (3.24) 

with some constant k. Here we first average over all possible motions of the vacancy 
up to time t {S(t)) for a fixed initial sea of dimers (with the vacancy in the center) and 
then average over uniformly chosen realizations of this initial dimer configuration. 

It is natural to compare this diffusion to that of a particle diffusing with the same 
dynamical rules on a uniformly chosen infinite spanning tree. This latter diffusion is 
characterized by a similar exponent r]Q with 

{S{t))o ~ kot^'^ , (3.25) 

where the average (■)o is now taken over all possible (infinite) spanning trees. 

A simple scaling argument relates the exponents 77 and 770. Denoting by s the size 
of the tree component of the spanning web for a fixed dimer configuration, we expect 
that S{t) ~ kof^^ with the infinite spanning tree diffusion exponent ryo as long as kof^" 
is less than s, while it saturates at S{t) ~ s for longer times. Using Eq. (|3.20|) for the 
distribution of sizes s, we can estimate 

{S{t)) dssp{s)+ ds kot"^' p{s) ^ t'^'>^^-^K (3.26) 

Jo ifcofo 
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Note that it is crucial that 1 < 5 < 2 for the two integrals to scale in the same way, 
determined by the t dependent limits of the integrals. In other words, the weak local- 
ization property of our model is essential for a well defined scaling. Taking 5 = 9/8, we 
have the relation 

v = lm, (3.27) 
which measures the lowering of the diffusion exponent due to weak localization. 



4. Monte Carlo simulations 
4-1. Simulations of spanning webs 

We have also found it revealing to perform numerical simulations of our model. As 
we shall see, we find very good agreement for all quantities for which we gave analytic 
predictions. Moreover, the simulations allow us to estimate both diffusion exponents i] 
and r]o independently, as well as other dynamical exponents. 

To begin with, we have simulated spanning webs using a modification of the Propp- 
Wilson cycle-popping algorithm for finding spanning trees of an arbitrary directed graph 
||18|| . This algorithm is presented in detail in Appendix B and works as follows: for every 
vertex of the odd white grid other than the vacancy vertex, we pick a possible outgoing 
arrow at random. If the resulting graph contains one or more cycles, we choose a cycle 
and "pop" it by picking a new random outgoing arrow for every vertex on the cycle. We 
then keep popping cycles until no cycle is left, resulting in a spanning tree. Propp and 
Wilson were able to show that this process terminates with probability one, and that 
the final result is independent of the order in which cycles are popped and produces 
trees uniformly distributed in the space of spanning trees. To generate spanning webs, 
we simply modified this algorithm by not popping cycles that encircle the vacancy. The 



proofs in Ref. |]T^ generalize straightforwardly to show that the modified algorithm 
results in web configurations uniform in the space of spanning webs. More precisely, 
this algorithm produces spanning webs with oriented loops, hence corresponding to a 
weight y = 2 per loop. It is then straightforward to correct for this weight so as to get 
an arbitrary value of y. 

As a preliminary check, we have confirmed that for a number of system sizes less 
than 50, the algorithm produces probabilities for the number of loops (in practice 0, 1 
or 2) in good agreement with the exact probabilities obtained for y = 2 by the methods 
of Sect. 3.1. 

We generated y — 2 spanning webs on odd white grids of size (2L-I-1) x (2L-I-1) with 
the vacancy at the center (these correspond to dimer grids of sizes (4L -|- 1) x (4L -|- 1)). 
Our system sizes were 2L+1 = 25, 35, 51, 75, 101, 151, 201, 301, 401, 501, 701, and 1001. 
The number of spanning webs generated was 25000 for all system sizes 2L -t- 1 < 301, 
30000 for 2L + 1 = 401, 10000 for 2L + 1 = 501, 3000 for 2L + 1 = 701, and 1500 for 
2L + 1 = 1101. The results were then transformed into probabilities for y = 4 spanning 
webs (dimer packings) by weighting configurations with C loops by an extra factor 2"^. 
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Fig. 16: Fraction of delocalized configurations, Vl, as a function of the 
system size. 



The delocahzation probabihty is measured by the fraction of spanning trees among 
spanning webs. This fraction is shown in Fig. 16 as a function of the system size and 
is weU fit by a power law. Dropping the two smallest system sizes (2L + 1 = 25 and 
2L + 1 = 35), which are still consistent with a power law fit but are too small to be safely 
in the large system limit, we get a delocahzation exponent 0.246 ± 0.006, in agreement 
with the value 1/4 predicted in Eq. ( p.l4|) . 

As for the distribution p{s), we first checked that the probability Pl(1) of a fully 
jammed state is in good agreement with both the exact value that we obtained for small 
system sizes [L < 50) and with the asymptotic prediction of Eq. ( |3.19| ) for p{l) (see 
Fig. 17). 

The distribution PL (s) for L = 200 is shown in Fig. 18. At this stage, it is important 
to analyze the effect of the boundary in a finite size grid, which can be described as 
follows: whenever the tree component of the spanning web touches the boundary, it must 
span the whole grid, hence has size (2L + 1)^. For s < L, the tree component cannot 
reach the boundary so we expect Pl{s) to be an accurate estimate of the asymptotic p{s). 
For L < s < (2L — 1)^, PLis) should become significantly lower than p{s) as s increases 
since configurations that contribute to Pl{s) are required to avoid the boundary, a 
constraint that does not exist in the asymptotic limit. Finally, the only possible value 
with s > (2L — 1)^ is s = (2L + 1)^ and corresponds to delocalized configurations, with 

a value pl (^(2L + l)^j = Vl that is not directly relevant to estimating p ^(2L + l)^j. 

This value is not shown in Fig. 18. By eye, the distribution appears to be a power 
law over almost four decades. Upon fitting the data however, a small but significant 
curvature appears for tree sizes greater than 10^'^, as expected from the above argument. 
To obtain the exponent, we thus fit only for tree sizes less than 200 (trees that cannot 
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Fig. 17 : Comparison between the exact finite size values of ^^(1), together 
with its conjectured exact asymptotic value, and the numerical probability 
of strict jamming. 
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Fig. 18: Distribution p{s) for the size s of the tree component as measured 
from spanning web simulations with free boundary conditions (green filled 
circles) on a 401 x 401 odd white grid (corresponding to a 801 x 801 dimer 
grid) with the vacancy at the center, and directly from dimer simulations 
with periodic boundary conditions (red open circles) on a 1101 x 1101 dimer 
grid. The fit (blue line) is for free boundary conditions and tree sizes less 
than 200. 
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reach the boundary), which gives a power law exponent 5 = 1.122 ± 0.008, in good 
agreement with the prediction 5 = 9/8 = 1.125 of Eq. (|3.23|) . It is surprising that 
the power law fit of the data extends all the way down to s = 1. In other words, the 
data are well approximated by a (normalized) pure power law distribution s-^/^/C(9/8). 
For s = 1, it gives a value 1/^(9/8) = 0.116. . ., to be compared with the exact value 
0.108... of Eq. (ICTI) . 

As we mentioned above, we expect that the asymptotic distribution ]?(s) should not 
depend on the finite size boundary conditions from which it is determined. In this vein, 
we have also generated dimer packings with periodic boundary conditions, for grids of 
odd linear size with a single vacancy. Note that the spanning web construction, central 
to this paper, no longer applies as, in particular, we have lost the global notion of distinct 
odd and even white grids in this case, or even the global notion of bicolorability. Still, 
these notions are preserved locally and for localized configurations with small enough s, 
the set of accessible sites again forms a tree whose statistics we can analyze along the 
same lines as for free boundary conditions. 

One advantage of periodic boundary conditions is that one can quickly generate 
many configurations, directly in the dimer setting, with the so-called "pivot algorithm" 
|]T^,^. Since there is no longer the constraint for localized configurations of avoiding 
the boundary, periodic boundary conditions have the further advantage of reaching 
larger values of s at fixed L, leading to smaller finite size effects. For a linear size of 
1101 (corresponding to L = 275 in the spanning web language), the resulting tree size 
distribution is well fit by a power law over a larger range (tree sizes 1 < s < 10^'^) 
than for spanning webs, with an exponent S = 1.121 ± 0.003, in good agreement with 
Eq. (ra . 

Looking at system sizes ranging from 101 to 1101, in steps of 100, we find that 
the fraction of delocalized configurations (now defined as configurations in which the 
vacancy can reach any vertex) has a power law exponent of 0.260 ± 0.005, again in 
agreement with Eq. (|3.14| ). 



4-2. Vacancy diffusion 

Let us now come to the prediction (|3.27|) relating the diffusion exponent r] for the 
growth with time of the number of sites visited by a diffusing vacancy to the correspond- 
ing exponent 770 for diffusion on an infinite spanning tree. To our knowledge, none of 
these exponents is known exactly and we therefore simulated both diffusion processes. 

For diffusion on uniform spanning trees, the simulations were done for free boundary 
conditions on a 401 x 401 grid, with the diffusing vacancy initially at the center of the 
lattice. At each time-step, the vacancy chooses one of the four compass directions at 
random and attempts a move in that direction. If there is a tree edge in that direction the 
move is carried out and otherwise, the vacancy stays at the same position. This slightly 
unusual random walk rule is chosen to mimic the standard Monte Carlo dynamics for 
the underlying dimers, where each time-step consists of an attempted dimer move. This 
means that a vacancy that has more possible moves should on average wait less time 
before moving. 
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Fig. 19: Number of distinct sites visited by a diffusing vacancy on a 
401 X 401 uniform spanning tree as a function of time. 
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Fig. 20: Squared displacement of a diffusing vacancy on a 401 x 401 
uniform spanning tree as a function of time. 

An ensemble of 1000 spanning trees was generated and, for each tree, the vacancy 
underwent 10^ time steps. The resulting graphs for the number of sites visited as a 
function of time, and the squared displacement as a function of time, are shown in 
Figs. 19 and 20. There is a relatively large curvature at smaller times [t < 100), so 
such times are excluded from the fit. Both graphs show good power law behavior by 
eye over the last five decades of time. Closer inspection shows that both graphs do have 
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a small but significant curvature, which we use to estimate the error bars of the slope 
(the statistical errors are negligible). 

From Fig. 19, we obtain r]o = 0.61±0.02, while from Fig. 20, we get Oq = 0.62±0.02, 
where we have defined the exponent 9o through 

(r2(t))ooct^° , (4.1) 

where r(t) is the Euclidean distance of the vacancy to the center of the grid at time t . 
These estimates are consistent with the two exponents being identical. 
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Fig. 21: Number of distinct sites visited by a diffusing vacancy on a 
401 X 401 uniform spanning web as a function of time. 

We then simulated the diffusion of a vacancy on spanning webs. The spanning 
webs were generated with the algorithm described in Sect. 4.1. An ensemble of 10000 
spanning webs was generated, each of the same size 401 x 401 as for the spanning trees 
above. For each web, the vacancy moved for the same total time and with the same 
dynamics as before. 

Figs. 21 and 22 show the number of sites visited and the squared displacement in 
the case of spanning webs. Again, we fit only for times t > 100, and use the curvature 
to estimate the error bars. Good power law behaviors are again seen over five decades. 
We obtain rj = 0.54 ± 0.03 and 9 = 0.56 ± 0.01 (with 9 defined as 9o but for spanning 
webs). Again, these estimates are consistent with the two exponents being identical. 

Finally, we estimate the ratio r]/r]o = 0.89 ± 0.06, in agreement with the prediction 
7/8 = 0.875 of Eq. (jO^ ). 

5. Discussion 

In this paper, we have analyzed the possible motion of an isolated vacancy in an 
otherwise fully packed dimer model on the square lattice. We find that the vacancy 
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Fig. 22: Squared displacement of a diffusing vacancy on a 401 x 401 
uniform spanning web as a function of time. 



exhibits weak localization: the size of its accessible domain scales as a power law with 
diverging average size. Exact finite size enumerations allowed us to identify several 
exponents as well as determine the probability for strict jamming. Some of these results 
are still awaiting exact and/or rigorous proofs. We also gave a universal relation between 
the exponents rj and 770 for the diffusion on spanning webs (vacancy diffusion) and that on 
spanning trees. We have no prediction for their individual values but from the numerical 
data, ryo is consistent with a simple value 5/8, which would result in a vacancy diffusion 
exponent rj equal to 35/64. 

The emergence of spanning webs as generalizations of spanning trees is quite natural 
here and one might hope that they will appear in other physical and mathematical 
problems. 

There are several directions in which the present work can be extended. Introduc- 
ing multiple vacancies will lead to vacancy interactions whose treatment will almost 
certainly require more elaborate geometrical structures. The power of the analysis in 
this paper relies heavily on special features of the square lattice. Extension to other 
lattices will require new insight. Key among them is the triangular lattice, for which 
we expect vacancies to be localized, with the size distribution of accessible sites now 
decaying exponentially at a rate determined by the known entropy mismatch between 
dimer configurations and spanning tree configurations. If we assume that this size is a 
good measure of the extent to which a vacancy perturbs its dimer background, the local- 
ization should be related to the exponential convergence of the "monomer-monomer" 
correlation function to a plateau at large separation (deconfinement) , as observed in 



Ref. [21|. Note finally that a random lattice version of the problem looks particularly 



promising for obtaining exact results. 
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Appendix A. Coulomb gas argument 



Let us present here a heuristic derivation of the relation ( p.l5| ) based on a Coulomb 
gas formulation of the problem. 
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Fig. 23: (a): The paths of the maze made of a spanning web (black) 
and its dual web (red) form a set of loops (green) that completely encode 
the (unrooted) spanning web. (b): The loops are fully-packed, concentric 
around the root of the original spanning web and can be oriented so as 
to follow a Manhattan orientation of the underlying grid (as indicated by 
arrows) . 

First we note that, as illustrated in Fig. 23, any spanning web configuration with 
C loop components can be alternatively coded by a configuration of 2£ -|- 1 loops on a 
square grid dual to the original dimer grid. These loops are such that: (i) each loop 
encloses the root vertex io, (ii) the loops are self- and mutually-avoiding and fully- 
packed, i.e. each vertex of the dual grid is visited by a loop, (iii) the loops can be 
oriented consistently so as to follow a Manhattan orientation on the dual grid. We 
can relax condition (i) and assign a weight uq for each loop that does not enclose io 
and a weight n for each loop that does enclose io. The correct statistics with a weight 
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y per loop component of the spanning web is then recovered by choosing no = and 
n = (and dividing by the residual weight associated with the loop around the tree 
component). The model is expected to lie in the universality class of the dense 0(no) 
model which can be described by a one-dimensional height field.0 In terms of this height 
field, the loop model is mapped onto a Coulomb gas model for which various exponents 
can be obtained exactly The anomalous weight n for the loops that enclose zq is 
properly accounted for by introducing an electric operator at with charge q such that 
n = 2 cos(7r(q' — uq)). The dimension of such an electric operator is 

X = , (A.l) 

where tiq = 2cos(7rwo) (0 < < 1), 6*0 = 1 ~ 'f^o and with the determination of 
q such that uq — 1/2 < q < uq + 1/2. Here we want tiq = 0, i.e. uq = 1/2 and 
n = = 2 cos{tiu) with u as in Eq. ( |3.15|) , i.e. q = u + 1/2. We end up with 

x = u' - -. (A.2) 
4 

The dimension x of the electric operator measures in particular the algebraic decay as 
of its average in a finite geometry of linear size L. For a square grid geometry, this 
average is nothing but the ratio Z{y;io)/Z{0;io) where the denominator is identified as 
the partition function of the system in the absence of electric operator. We immediately 
deduce the relation x = 7(2/) — 7(0), i.e. 

7(2/) = «'-^ +7(0) . (A.3) 

To end the argument and recover Eq. (|3.15| ), we simply rely on a direct calculation of 
7(0) = 1/2, as obtained in Ref. [TO from the exact product formula for Z{0; zq) = ^tree- 



Appendix B. Generating spanning webs via cycle-popping 

In this appendix we explain how to generate a random spanning web using a "cycle- 
popping" algorithm inspired from the Propp- Wilson algorithm for the generation of a 

^ It is known that fully-packed loops are in general in a different universality class than 
dense loops. In particular, on the square lattice, the fully-packed loop universality class is 
described by a 3-dimensional height variable. For the model at hand with the extra constraint 
(iii) above, however, two of the height components are eliminated, leading to an effectively 
one-dimensional model. 
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(a) 



(b) 



Fig. 24: An example of a configuration (a) obtained by drawing at random 
an edge from each vertex but the root io. A possible configuration (b) 
obtained from (a) by popping the magenta cycle. The tree component 
fiowing to io (here in green) can only grow in the process. To generate 
spanning webs instead of spanning trees, we make all cycles winding around 
io (such as the brown cycle) unpoppable. 

random spanning tree. More precisely we extend to random webs the RandomTreeWith- 
RootO procedure explained in Sect. 6 of Ref. [ |TB| . 

Let us begin by recalling how the Propp- Wilson algorithm works in our specific 
setting. Start with a rectangular grid on which the spanning tree is to be constructed, 
with the root io at a given position. At each vertex distinct from io, draw a random 
incident edge uniformly and independently from the other vertices. Pictorially the edge 
is marked and oriented away from the vertex : a possible outcome of this procedure is 
illustrated in Fig. 24-(a). Note the similarity with the proof of the determinant formula 
(^]2|). By chance, the graph made of the selected edges can be a spanning tree rooted 
at io (with edges oriented towards the root), but more likely it will contain one or 
more (disjoint) cycles (loops) as in Fig. 24-(a). Cycle-popping consists in choosing an 
arbitrary cycle and popping it, i.e. for each vertex on the cycle, drawing a new random 
outgoing edge uniformly and independently of all previous draws. For instance. Fig. 24- 
(b) shows a possible outcome after the popping of the magenta cycle in Fig. 24-(a). 
The resulting graph can again contain cycles (some preexisting and some created by 
popping): in that case repeat the procedure by choosing another cycle and popping 
it. Otherwise the resulting graph is a spanning tree and the procedure terminates. In 
Fig. 24, the green edges "fiow" to the root, hence they cannot be in a cycle and they 
must belong to the final tree. One then sees that the green graph can only grow at each 
step, eventually covering the whole grid. 

There are several ways to choose which cycle to pop at each step of the procedure 
(one could imagine doing it deterministically or randomly), but Propp and Wilson have 
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shown that the actual procedure is (in some precise sense) irrelevant. Furthermore: (i) 
the procedure terminates almost surely, i.e. only a finite number of cycles have to be 
popped before none remains, and (ii) the resulting spanning tree is a perfect sample 
of the uniform measure on the set of spanning trees rooted at io. Propp and Wilson 
actually deal with a slightly more general case, not needed here, and compute bounds 
on the running time, which establish that this is a very efficient algorithm. 

The modification of the algorithm to generate spanning webs is easy. Perform the 
procedure as before, except that cycles winding around the root io are now considered 
unpoppable and stay in the resulting graph. For instance in Fig. 24- (b) the brown cycle 
would be unpoppable if we were to generate a spanning web. Only the cycles not 
winding around the root are popped, until none remains. Propp and Wilson's analysis 
(in Sect. 7 of Ref. [0) extends straightforwardly to the case of the unpoppable cycles, 
and we find that : (i) the procedure terminates almost surely, and (ii) the resulting 
graph is a now a perfect sample of the uniform measure on graphs (made out of oriented 
edges, with exactly one edge going out of every vertex, except for the root that has 
none) whose cycles are all unpoppable. One sees easily that such a graph is nothing 
but a spanning web rooted at io, as defined in Sect. 2.2, with additional orientations 
such that each edge not belonging to a cycle "fiows" towards the root or an attractor 
cycle. Since each cycle has two possible orientations that are equally likely, we conclude 
that the probability to obtain a given unoriented spanning web containing C loops is 
2'^/Z(2,zo), where Z(2,io) is the generating function for spanning webs on the grid at 
hand, rooted at io and counted with a weight y = 2 per loop, as defined in Eq. (|2.4| ) . This 
is to be contrasted with the weight y = 4 per loop that arises from the correspondence 
with dimer configurations with a vacancy, drawn with uniform probability. Knowing 
the exact bias, it is however straightforward to translate the statistical properties of 
spanning webs as measured through this algorithm into dimer statistics. Moreover, 
one can correct the bias directly in the algorithm, at the price of introducing suitable 
weights for edges along a seam. 

In conclusion, we have provided an algorithm for the generation of random spanning 
webs belonging to the class of exact or perfect algorithms [p!8| , p3[| . By the correspondence 
of Sect. 2.2, it can be used to simulate dimer configurations on a square lattice with 
a vacancy. We have not performed a detailed analysis of its efficiency, but we believe 
it is comparable to the original Propp- Wilson algorithm for rooted spanning trees, and 
much better than generic CFTP-type algorithms. The comparison with other Monte 



Carlo algorithms, such as the pivot algorithm [ p!9| , p0| for the case of periodic boundary 
conditions is not so clear: typically, these are designed with efficiency in mind, but 
without an exact knowledge of their "randomness". Both approaches provide useful 
and mutually consistent results in our study. 
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